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ABSTRACT 

We use the billion-particle Hubble Volume simulations to make statistical predictions 
for the distribution of galaxy clusters that will be observed by the Planck Surveyor 
satellite through their effect on the cosmic microwave background — the Sunyaev- 
Zel'dovich effect. We utilize the lightcone datasets for both critical density (tCDM) 
and flat low-density (ACDM) cosmologies: a 'full-sky' survey out to z ~ 0.5, two 
'octant' datasets out to beyond z = 1 and a fOO square degree dataset extending to 
z ~ 4. Making simple, but robust, assumptions regarding both the thermodynamic 
state of the gas and the detection of objects against an unresolved background, we 
present the expected number of SZ sources as a function of redshift and angular size, 
and also by flux (for both the thermal and kinetic effects) for 3 of the relevant HFI 
frequency channels. We confirm the expectation that Planck will detect around 5 x 10 4 
clusters, though the exact number is sensitive to the choice of several parameters 
including the baryon fraction, and also to the cluster density profile, so that either 
cosmology may predict more clusters. We also find that the majority of detected 
sources should be at z < 1 .5, and we estimate that around one per cent of clusters 
will be spatially resolved by Planck, though this has a large uncertainty. 

Key words: methods: numerical - cosmology: theory 



One of the most exciting techniques for surveying for galaxy 
clusters at high redshift, which are a powerful probe of cos- 
mology, is via the Sunyaev-Zel'dovich (SZ) effect (Sunyaev 
& Zel'dovich 1972, 1980; Rephaeli 1995; Birkinshaw 1999). 
The effect is readily measured, and indeed mapped, in known 
clusters (e.g. Carlstrom et al. 2000) and several experiments 
are in development which will be able to make surveys. The 
two main survey types are ground-based experiments prob- 
ing small fractions of the sky at high angular resolution, 
including interferometers such as BIMA and AmiBA and 
single-dish experiments such as Viper and BOLOCAM, and 
satellite projects making all-sky surveys at lower resolutions, 
such as the Planck Surveyor satellite (hereafter referred to 
as Planck). 

The SZ effect is the distortion in the spectrum of the 
cosmic microwave background (CMB) as a result of the 
change in energy of the CMB photons when they scatter 
from hot gas, usually that in galaxy clusters. The effect is 
typically separated into two components: the kinetic SZ ef- 
fect is the Doppler shift from the bulk motion of the clus- 
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ter, while the thermal SZ effect is the gain in energy that 
the photons experience because the gas is at a much higher 
temperature than the CMB. As well as being of interest in 
its own right, a detailed understanding of the SZ effect is 
necessary in order to remove its effect on primary cosmic 
microwave background anisotropies. 

To date, two main approaches have been taken to model 
the SZ effect: semi-analytical and direct numerical sim- 
ulation. The former method ordinarily assumes that the 
mass distribution of haloes at a given redshift is accurately 
represented using the Press-Schechter formalism (Press & 
Schechter 1974). Making simplifying assumptions such as 
isothermality allows the abundance of haloes as a function of 
the thermal SZ signal to be constructed (e.g. Barbosa et al. 
1996; Eke, Cole & Frenk 1996). Although Press-Schechter is 
in good agreement with results from dark matter simulations 
over a wide range of mass scales, it fails for the rare, most 
massive, objects in a way which has now been carefully quan- 
tified (Jenkins et al. 2000), and indeed these corrections have 
been incorporated in a recent study by Bartelmann (2000). 
Such methods have the drawback of making simple assump- 
tions about the structure of the gas within the haloes, and 
also cannot be used to study the kinetic effect since they 
do not allow a direct prediction of the halo peculiar veloc- 
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ity field. Although it has recently become possible to study 
the kinetic effect through more detailed semi-analytic tech- 
niques (e.g. Benson et al. 2000; Valageas, Balbi & Silk 2000), 
such treatments are extremely complex. 

Numerical simulations, though much more computa- 
tionally intensive, are able to directly model the formation of 
structure in both baryons and dark matter. Hydrodynamical 
iV-body simulations are the natural tool; early studies were 
made by Thomas & Carlberg (1989) and by Scaramella, Cen 
& Ostriker (1993), and more recently the much-improved 
numerical technology has been exploited by several groups 
(da Silva et al. 2000a,b; Refregier et al. 2000; Seljak, Bur- 
well & Pen 2001; Springel, White & Hernquist 2000; Gnedin 
& Jaffe 2000). Such work is now reaching a state of matu- 
rity where the importance of additional physics is starting 
to be considered. These simulations have proven capable of 
making maps of the SZ effect at high resolution (around one 
arcminute in most of the quoted papers, but at arcsecond 
resolution in the small-scale simulations of Gnedin & Jaffe 
2000), but it has not been possible to make maps covering 
a large area of sky and hence capturing the rarest, brightest 
objects. 

In this paper we take an approach which can be consid- 
ered a mid-point between Press-Schechter and full hydro- 
dynamical simulations, which is to use the largest existing 
iV-body simulations. These are the Hubble Volume simula- 
tions, run by the Virgo Consortium (Evrard et al. 2001, in 
preparation; Jenkins et al. 2000; Colberg et al. 2000), for two 
popular cosmological models — a critical-density model with 
a modified power spectrum (which we refer to as rCDM) and 
a flat, low-density model with a non-zero cosmological con- 
stant (which we refer to as ACDM) . The main advantage of 
these simulations over hydrodynamical simulations is their 
volume; typical hydrodynamical simulations are performed 
using a box-length of ~ 100 h~ x Mpc; the Hubble Volume 
simulations are over an order of magnitude larger; both mod- 
els contain N — 10 9 particles, within a 2/i _1 Gpc (comov- 
ing) box for the rCDM simulation and a 3/i _1 Gpc box 
for the ACDM simulation. In particular, this allowed "light- 
cones" to be constructed as the simulations were performed, 
obviating the need to stack simulation boxes as in the hy- 
drodynamical simulations approach. Hence, the provision of 
lightcones allows us to measure the abundance of clusters 
over the range of appropriate redshifts, for a range of mass 
scales that is larger than in any of the previous work. Work- 
ing with the simulation data directly avoids the need for 
approximate mass functions such as Press-Schechter, and 
also gives us direct access to the kinetic effect. However, 
the fact that the simulations are dark matter only means 
that we have to make some of the simplifying assumptions 
concerning cluster properties characteristic of semi-analytic 
methods. Therefore, such large-scale simulations are appro- 
priate for making predictions for all-sky surveys rather than 
high-resolution measurements, and so our main focus of this 
paper is directed towards predictions for the Planck satel- 
lite. As well as source counts, we will be interested in the 
expected size and redshift distributions of the clusters. 



2.1 SZ definitions 

The SZ effect produces a fluctuation in the surface bright- 
ness of the CMB from the inverse Compton scattering of the 
CMB photons off electrons within an ionized plasma. This 
is due to the internal (thermal) motions of the electrons 
and also the bulk (kinetic) motion of the plasma, relative 
to the CMB reference frame. The total fluctuation can be 
expressed as (e.g. Rephaeli 1995) 



SIu = Io [g(x)y - h(x)Pr] 



(1) 



where Jo = 2fc|To / h 2 c 2 , To = 2.725 K is the mean tempera- 
ture of the microwave background (Mather et al. 1999) and 
(5 = v r /c is the radial component (i.e. projected along the 
line of sight) of the plasma's bulk velocity in units of the 
speed of light in a vacuum. The functions g(x) and h(x), 
known as spectral functions, contain frequency-dependent 
information for the thermal and kinetic effect respectively, 
with x — hu/kBTo- Specifically, 



and 



e x + 1 



h(x) = 



(e* - l) 2 



(2) 



(3) 



The dependence of 5I V on the thermodynamic state of the 
gas is given by the Comptonization parameter, y, and the 
optical depth due to electron scattering, r. For the thermal 
effect 



y = / n e l e dl , 



(4) 



i.e. an integral of the electron pressure along the line of sight, 
while for the kinetic SZ effect we have 



r = o"t 



n e dl . 



(5) 



The total SZ fluxQ of a source is then the integral of 8I U 
over its subtended solid angle 



S v {x) = So / [g{x)y- h(x)/3r] dft, 



(6) 



where So 



2.29 x 10 4 mJy arcmm 



(lmjy 



10 _26 ergs _1 m _2 Hz _1 ), assuming that the dimensions of dQ 
are arcmin 2 (the bracketed terms are dimensionless) . Note 
that y and r are positive quantities, whereas (3 can be pos- 
itive (if the source is moving towards the observer) or neg- 
ative (if the source is moving away from the observer). The 
function g(x) is zero for x ~ 3.83 (y ~ 217 GHz), nega- 
tive if x < 3.83 causing a decrement in the CMB surface 
brightness, and positive if x > 3.83, causing an increment. 
The function h(x) is always positive. Hence, the signal due 
to the SZ effect can be deconvolved from the total CMB 
signal, into thermal and kinetic components if the signal is 
observed over a range of frequencies. 

For the thermal effect, it is convenient to quote the 
frequency-independent part of equation ([]) 



2 METHOD ^ e use * ne term "flux" to mean the change in flux caused by 

the SZ effect relative to the mean flux of the CMB. 
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Table 1. Key parameters of the Hubble Volume simulations. 



Model 


fi 


n A 


h 


r 


0"8 


rCDM 
ACDM 


1.0 
0.3 


0.0 
0.7 


0.5 
0.7 


0.21 
0.21 


0.6 
0.9 



m/10 12 ft _1 M© 



Y = d ~A j ydA, (7) 

where d£l = dA/d A , cLa is the (cosmology-dependent) an- 
gular diameter distance to the source and the integral is 
performed over the source's projected area. Since y is dimen- 
sionless, the dimension of Y is that of a solid angle. Hence, Y 
represents a weighted measure of the solid angle of the CMB 
that is obscured by the given source. Note that although y 
is independent of cosmology (except through the expected 
distribution of sources) , Y depends on cosmological param- 
eters via the relation between the (inferred) physical and 
(observed) angular size of the source, through cLa- 

2.2 Simulation details and object selection 

We consider both cosmological models adopted for the Hub- 
ble Volume simulations, labelled rCDM and ACDM re- 
spectively. Key parameters of the simulations are listed 
in Table [l], namely the density parameter Qo, the cosmo- 
logical constant Qa = A/3Hq, the Hubble constant h = 
Ho/ W0 kms -1 Mpc -1 , the power spectrum shape parame- 
ter F, the variance as of linear fluctuations smoothed with 
a top-hat filter on the (comoving) scale of 8 h~ x Mpc, ex- 
trapolated to z = 0, and the particle mass m. 

We use the publicly- available dark matter halo cata- 
logues from the lightcone datasets.f[] Haloes are denned as 
overdense spheres of particles with a mean internal density, 
(p) = A c p CT (z), where p C r(z) = 3H(z) 2 /8ttG is the crit- 
ical density and H(z) is the Hubble parameter. For these 
catalogues A c = 200, which for Q — 1 is close to the canon- 
ical value of 178 predicted by the spherical top-hat collapse 
model (Lacey & Cole 1993). However, for the ACDM model 
our adopted value of 200 is around a factor of 2 larger than 
predicted by the top-hat model at z = 0. For the main re- 
sults presented in this paper, the choice of density threshold 
is only used to normalize the mass distribution within haloes 
and not to explicitly define their sizes. 

The information extracted from the simulations is the 
mass, redshift and velocity distributions of the dark mat- 
ter haloes. The range of mass scales for the haloes resolved 
by the Hubble Volume simulations is nearly two orders of 
mag nitude (~ 3 x 10 13 to 10 15 /i _1 M Q ), covering the scales 
relevant to galaxy clusters. We assume that the brightest SZ 
sources are due to the intracluster medium (ICM) trapped 
within the potential wells of cluster haloes, known to be a ro- 
bust approximation from the results of full hydrodynamical 
simulations. 

We use the lightcone datasets specific to each cosmology 
to construct two mock cluster surveys. Firstly, we utilize 
the MS sphere dataset, covering the full 4n steradians from 
z = to z = 0.44 for rCDM and z = 0.58 for ACDM. 

t These are available at www.mpa-garching.mpg.de/Virgo/ 
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Table 2. Key properties of the lightcone datasets. 



Model 


Dataset 


2-min 


ax 


fi/sr 




TLvJJJyl 


Mo 


U.UU 


n A A 
U.44 


1 — 

171 


llOOOOO 


rCDM 


NO 


0.44 


1.25 


tt/2 


227525 


rCDM 


PO 


0.44 


1.25 


tt/2 


225820 


tCDM 


DW 


1.25 


4.37 


0.03 


199 


ACDM 


MS 


0.00 


0.58 


4-7T 


1539281 


ACDM 


NO 


0.58 


1.46 


tt/2 


625138 


ACDM 


PO 


0.58 


1.46 


tt/2 


612549 


ACDM 


DW 


1.46 


4.60 


0.03 


2689 



Beyond these redshift limits, we graft on both PO & NO 
octant datasets, spanning a total of ir steradians between 
2 = 0.44 (z = 0.58) and z = 1.25 (z = 1.44) for tCDM 
(ACDM) |U Finally, we use the DW datasets which cover 
100 square degrees and extend to z — 4.37 for rCDM and 
z — 4.6 for ACDM. Combining the datasets allows us to 
have a mass-limited survey beyond the redshifts where the 
first resolved structures in the simulations form. 

We summarize key properties of each dataset used in 
Table |^. Columns 1 & 2 give the cosmological model and 
dataset label respectively; Columns 3 & 4 list the redshift 
limits used in the surveys; Column 5 gives the solid angle 
of each dataset and Column 6 lists the number of objects 
extracted from each dataset. 

2.3 Baryon distribution 

Given the mass and redshift of each halo, calculating the 
resulting SZ fluxes requires knowledge about the thermody- 
namic (both density and temperature) state of intracluster 
gas. Although the dimensions of the simulations are ideal 
for our purposes, they did not explicitly model a baryonic 
component. Therefore, we need to make some of the same 
assumptions made in analytical calculations. The first of 
these is to define the relationship between the electrons and 
baryons, and for this, we assume that the ICM is fully ion- 
ized, so that n e = 0.88 p/mn (assuming a helium mass frac- 
tion, Y — 0.24), where ?tih is the mass of a hydrogen atom 
and p is the baryon mass density. 

To relate the baryonic mass to the total mass of each 
cluster, we assume a global baryon fraction, / b = M b /M to t- 
We neglect the presence of cold gas and stars, which only 
reduces the gas fraction by around 20 per cent for the mass 
range of clusters relevant to this study (e.g. Balogh et al. 
2001). The SZ signal simply scales proportional to fh, so 
theoretical predictions can readily be rescaled to a differ- 
ent value (though ultimately there are nontrivial effects in 
deciding whether extended sources are detectable against a 
background). As emphasized recently by Bartelmann (2000), 
the choice of baryon fraction can make a significant differ- 
ence when contrasting cosmologies. The two main observa- 
tional constraints of /b are from X-ray observations of clus- 
ters (e.g. White et al. 1993; Ettori & Fabian 1999) and from 

§ We note that there is a small degree of overlap between the PO 
& NO datasets at the highest redshifts. However, removing one 
of these datasets does not significantly affect our results and so 
we include both for the sake of completeness. 
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primordial nucleosynthesis calculations (e.g. Buries & Tytler 
1998), assuming fb = Slb/^o and that the mix in clus- 
ters represents the cosmic mean. Both determinations are 
concordant with a low-density cosmology but discrepant if 
there is a critical density (the so-called baryon catastrophe) . 
Therefore, it appears that if critical density is to have any 
chance of fitting the whole range of available observations, 
then the baryon density must be higher than nucleosynthesis 
allows, and we note in passing that there are some moder- 
ately well motivated modifications of standard nucleosynthe- 
sis which make the element abundances compatible with a 
higher baryon density (see e.g. Kang & Steigman 1994; Les- 
gourgues & Peloso 2000; Kaplinghat & Turner 2001). We 
therefore adopt the observed baryon fraction from clusters 
when plotting results, taking the value from Ettori & Fabian 
(1999) of fb = 0.06 h- 3/2 . This gives f b = 0.10 for the low- 
density case and fb = 0.17 for the high-density case. The 
former is in reasonable agreement with the preferred stan- 
dard nucleosynthesis value fib = 0.019/i~ 2 (Buries & Tytler 
1998). In the critical-density case our main results will as- 
sume the cluster baryon fraction number, but occasionally 
we will show results based instead on the nucleosynthesis 
value fb = 0.076 for comparison. 



2.4 Cluster temperatures 

To evaluate the integral in equation (^) requires knowledge 
of the electron pressure profile of the cluster. To simplify 
the discussion, we assume that the ICM can be modelled 
as isothermal for all but the outer regions of clusters. T e 
can then be removed from the integral in equation (^), and 
in equation (^) the thermal component scales proportional 
to the total energy of the cluster and the kinetic as the 
cluster mass. Hydrodynamical simulations have confirmed 
that radial variations in the ICM temperature are weak, 
being within a factor of two out to the virial radius (e.g. 
Eke, Navarro & Frenk 1998; Pearce et al. 2000). We have 
calibrated temperature to mass using the relation 



M 200 = 2.5 x 10 13 (Tx/keV) 1 ' 5 (1 + z)" z/2 h' 1 M Q . 



(8) 



where M200 is the mass within each halo where the mean in- 
ternal density, (p) = 200p cr , as defined above. Equation (^) 
is consistent with the observed relation derived by Horner, 
Mushotzky & Scharf (1999), using spatially-resolved ASCA 
X-ray temperature data for a cluster sample with z < 0.1. 
(Note that we assume equivalence between the observed 
emission-weighted temperatures and the intrinsic electron 
temperature of the gas.) To account for non-zero redshift 
scaling, we assume that T oc M200 /R200 (the virial relation) 
and 



R 200 = I |g M 200 



200 H 2 {z) 



1/3 



(9) 



We note that H~ 2/3 (z) oc Q} ,z (z)/{1 + z) and neglect the 
weak dependence on fi. This results in the presence of the 
(1 + z) -3 / 2 term on the right-hand side of equation (H). 



2.5 Density profiles 

For our main results, we consider two forms for the distri- 
bution of baryons within cluster potential wells, since the 



precise form of p(r) is not yet well determined from X- 
ray observations (although significant improvements should 
come from Chandra and XMM-Newton observations). His- 
torically, cluster ICM density profiles are obtained by fitting 
a P model to the (azimuthally-averaged) observed surface 
brightness distribution and making an assumption regarding 
the temperature distribution (usually isothermality) . How- 
ever, most of these observations suffer from the inability to 
probe the outer radii of clusters because of the high X-ray 
background (which is a combination of instrumental noise 
and the intrinsic X-ray background of unresolved sources), 
particularly in low temperature systems. Consequently, the 
outer slope of the density profile is usually ill-constrained. In 
light of these difficulties, we choose to consider two density 
profiles for the gas that are already well documented in the 
literature. The first of these is the /3 model profile, the stan- 
dard model for fitting the observed ICM surface brightness 
distributions outwith cooling flow regions 



Pp{ r ) 



ps(.o) 



[l + (r/r c ) 2 ]-W 2 ' 



(10) 



where the parameters r c and f3 determine the core radius and 
the asymptotic slope of log(p) at large radii respectively. We 
calculate P/3(0) assuming that the integrated baryonic mass 
of the cluster should equal fb M200 at 7?2oo- For the core 
radius, we set r c = 0.07i?2oo (e.g. Lloyd-Davies, Ponman 
& Cannon 2000). Measured f3 values for individual clus- 
ters are typically between 0.5 and 0.8 (e.g. Jones & Forman 
1984,1999; Horner et al. 1999; Mohr, Mathiesen & Evrard 
1999; see also Vikhlinin, Forman & Jones 1999). For the re- 
sults presented in this paper we adopt (3 = 5/6 (~ 0.83) in- 
ferring an asymptotic slope at large radii of —2.5. Although 
this is rather high, we note that our assumption of isother- 
mality is problematic when calculating y for shallow density 
profiles, since the temperature should drop significantly out- 
side the virial radius (this is verified by simulations) . There- 
fore, we choose to interpret our choice of (5 as a way in which 
to parameterize the radial fall-off of the ICM pressure, p T, 
at large radii. 

The second density profile considered in this paper is 
well-motivated by results from numerical simulations: the 
NFW profile (Navarro, Frenk & White 1997), which provides 
a good fit to average dark matter halo density profiles over 
a wide range of radii, for all popular CDM and scale-free 
cosmological models 



PNFw(r) 



Pnpw(0) 



(11) 



The NFW model contains just one adjustable parameter, 
known as the concentration parameter, c = R2oo/r s and 
determines the radius at which the density profile steepens 
from a slope of —1 to an outer slope of —3. We normalize 
the profile as for the (5 model and calculate c as a function 
of halo mass and redshift using the method described in 
Appendix A of Navarro et al. (1997). 

To summarize, the inner regions of clusters are mod- 
elled as isothermal with either a j3 model (constant core) or 
NFW (cuspy) density profile. In the outer regions of clus- 
ters, the NFW profile is sufficiently steep that variations in 
temperature are unimportant, whereas a value of (3 has been 
chosen in order to represent what we anticipate as the weak- 
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est dependence of y with radius. Hence, we argue that the 
two models should span a reasonable distribution of what 
the true SZ surface brightness profiles are likely to be. 

2.6 SZ Fluxes and object sizes 

We concentrate on calculating fluxes for 3 frequencies: 143, 
217 and 353 GHz. These frequencies correspond to the cen- 
tral values of 3 channels on the Planck HFI detector (Puget 
et al. 1998) . Although there are other channels at higher fre- 
quencies (545 and 857 GHz), the signal at those frequencies 
will be dominated by dust emission. These frequencies have 
been specifically chosen in order to allow the best separa- 
tion of the signal due to the SZ effect and the true primordial 
temperature fluctuations. 143 GHz corresponds to the max- 
imum decrement in the CMB spectrum that is induced by 
the thermal SZ effect. At 217 GHz, there is no signal due 
to the thermal SZ effect, allowing the contribution from the 
kinetic effect to be measured. 353 GHz corresponds to the 
maximum increment in the CMB spectrum induced by the 
thermal SZ effect. 

The Planck HFI detector will not measure tempera- 
ture fluctuations only at these precise frequencies. Rather, 
each channel will have its own specific frequency response 
function. We model each response function as a Gaussian 
and assume that the full- width half-maximum (FWHM), 
1(y — Vq)/vo = 0.37 (Hobson et al. 1998), where vo is the 
central frequency. We calculate the average fluxes that are 
measured for each central frequency by replacing g(x) (equa- 
tion ^ and h(x) (equation |^) with 

G(x ) = -±- J g(x) exp (- {x "fQ dx, (12) 
and 

H(x ) = -J=- J h{x) exp (- ( X ~J 0) J dx, (13) 

where a = 0.37a;o/2\/ln 2. For an ideal detector of width 
v— vo — > 0, (i.e. measuring the flux at a single frequency), the 
response function would approach the Dirac delta function 
and g(x) and h(x) would both be recovered. 

The simplest model used to construct SZ fluxes is to 
assume that the size of each halo is -R200, in which case the 
total flux (equation |^) does not depend on the thermal dis- 
tribution of the intracluster gas. For the thermal SZ effect, 
we can define the integrated ly out to -R200, Y200, by combin- 
ing equations (JiJ), (0) and (^) to obtain 

where M200 is in units of /i -1 Mq and dA is in units of 
h^ 1 Mpc. This model assumes that there is negligible sig- 
nal from haloes outside -R200 , which would be the case if the 
gas is only ionized inside R200 (i.e. with a temperature sig- 
nificantly greater than 10 4 K). However, this approach does 
not allow us to account for beam convolution, which requires 
the exact distribution of the gas to be modelled, and so we 
choose to present this model only to contrast with our main 
results. 

The total SZ flux from each cluster is calculated by 
first defining a beam-convolved y profile, assuming that the 
Planck beam profile is adequately represented by a Gaussian 



Vc (0)= / y{e')b{\e-e'\)d 2 e , (15) 

where 9 = \0\ is the projected angular radius of the cluster 
and 

b(\0-e'\) = -Lexp(- 10 'J |2 V (16) 

Standard deviations are calculated using the full-width half- 
maxima quoted in the literature, i.e. approximately 8 ar- 
cminutes for the 143 GHz channel and 5 arcminutes for the 
217 and 353 GHz channels. The total SZ flux is then 

Sc{x ) = So / 2nd [G{x )y c (6) - H(x )Pr c (9)] d6, (17) 
Jo 

where the beam-convolved optical depth, r c oc y c /T in the 
isothermal limit. The upper limit in equation (|17|), 9h s is de- 
fined as the angular radius where y c equals the background 
due to unresolved sources, y^ g . Values of j/b g are not easily 
established and one must account for not only the mean y, 
(y) but also fluctuations about (y) , which will depend on the 
beam-size. At present, only an observed upper limit exists 
to (y), of (y) < 1.5 x 10~ 5 from the FIRA S instrument on 
the COBE satellite (Fixsen et al. 1996), while fluctuations 
about (y) are constrained only on very large angular scales. 
However the background can be predicted from both the 
Press-Schechter approach (Bartelmann 2000) and from hy- 
drodynamical simulations (da Silva et al. 2000b) . We adopt 
values from the latter and choose to define the background, 
2/bg as the 3 sigma fluctuation above (y). This value is ap- 
proximately ybg = 5 x 10~ 6 for a Gaussian beam with a 
FWHM of 5 arcminutes. 

Finally, to estimate the sizes of objects that will be ob- 
served by Planck, we cannot simply take the #b g values, 
since this does not give information regarding whether a 
source is either a point source (which may be extended due 
to the beam convolution) or actually resolved resolved by 
the beam. Instead, we define the angular size of each object, 
9 S , to be the angular radius at which y c (9) = y c (0)/2, i.e. 
the half-width half-maximum. For objects where 9 S > 6>b g 
(i.e. j/bg > j/ c (0)/2), then 9 3 is set to #b g - We define sources 
with 9 S greater than the beam FWHM as resolved, since this 
will produce flux in more than one detector pixel (which will 
typically have a size of order the beam FWHM). 



3 RESULTS 

3.1 Halo distributions 

Fig. jj] shows the number of haloes per unit redshift per 
square degree, as a function of redshift, for all haloes ex- 
tracted from both simulations. Firstly, the area under the 
ACDM redshift distribution is greater than for rCDM due 
to the larger number of objects in the ACDM catalogue 
(~ 2.8 x 10 6 , as opposed to ~ 1.6 x 10 6 ). Secondly, the peak 
of the ACDM distribution is higher than for the rCDM case. 
For example, the rCDM simulation predicts that the num- 
ber of objects per unit redshift is around 100 per square de- 
gree, whereas the ACDM simulation predicts that the peak 
number is around 50 per cent larger. Secondly, the peak 
of the rCDM distribution occurs at z ~ 0.4, whereas the 
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Figure 1. Rcdshift distributions of haloes from each survey. 
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Figure 2. Halo mass functions for the rCDM survey (top panel) 
and ACDM survey (bottom panel). The solid line is for the com- 
bined datasets and the dotted/dashed lines illustrate contribu- 
tions from individual datasets covering various redshift intervals. 

peak of the ACDM distribution occurs at z ~ 0.8, result- 
ing in a much broader redshift distribution. Consequently, 
the tCDM survey extends to z ~ 2 whereas the ACDM 
survey extends to z ~ 3. These differences arise due to the 
slower rate of structure formation in the low-density case, 
which requires a larger normalization of density fluctuations 
(through as) in order to produce the same z — abundance 
of rich clusters as in the critical-density model. 

Differential mass functions of haloes are illustrated in 
Fig. I for rCDM (top panel) and for ACDM (bottom panel). 



The solid line in each figure illustrates the mass func- 
tion for the complete survey, in addition to contributions 
from the MS dataset (dotted line), the combined PO & 
NO dataset (short-dashed line) and the DW dataset (long- 
dashed line). Hence, each component illustrates the contri- 
bution from 3 redshift intervals. The majority of high-mass 
objects (> 1O 15 /i _1 M ) are at low redshift and the majority 
of haloes at lower masses come from intermediate redshift. 
This is due to the characteristic halo mass (known as M* , 
whose value roughly coincides with the knee of the halo mass 
function) decreasing with increasing redshift, as is the case 
in all hierarchical models. 

We also note that the haloes with redshifts above z ~ 
1.2 — 1.5 (i.e. from the DW dataset) make only a small con- 
tribution to the total mass function in either model. The 
amplitude of the DW mass function is at least a factor of 5 
lower than the total mass function in the ACDM case, and 
over a factor of 20 lower in the tCDM case, where structure 
evolves more rapidly with redshift. 

The cumulative number of haloes above Y200 is plot- 
ted in Fig. ^| We note that our datasets are mass-limited 
and not flux-limited, as is the case with real surveys, and 
so if we choose to consider the full mass range of haloes, 
we may be incomplete in flux due to the absence of smaller 
mass sources that are sufficiently close to z = to have 
Y values comparable to our resolved sources. As noted by 
Bartelmann (2000), a conservative sensitivity limit of forth- 
coming Planck surveys is around Yn m = 3 x 10~ 4 arcmin 2 . By 
considering only haloes with Y200 > Ylm, we find that 90 per 
cent of these objects have masses greater than 25 particles 
for tCDM and 40 particles for ACDM respectively (the mass 
limit is equivalent to 12 particles). Hence, our constructed 
surveys should be sufficient to describe the distribution of 
objects with Y > 3 x 10 _4 arcmin 2 . The vertical dashed line 
in Fig. ^ denotes Yu m . 

The solid and dashed lines in Fig. ^ represent the source 
counts above a given Y200 for the rCDM and ACDM surveys 
respectively, assuming a value of /b calculated from clusters. 
The number of sources declines more rapidly for the ACDM 
case than for the rCDM case. For Y200 > 10 _4 arcmin 2 both 
surveys predict approximately 10 clusters per square degree, 
but for Y200 > O.larcmin 2 the rCDM survey predicts around 
an order of magnitude more objects than the ACDM survey. 
This is because the most massive clusters in the ACDM 
model are, on average, at higher redshift than those in the 
rCDM model. Furthermore, the angular diameter distance 
at fixed redshift is larger for the ACDM model than for the 
rCDM model. These two differences conspire to produce 
smaller Y200 values for the largest clusters in the ACDM 
model than in the rCDM model. However, for lower Y200 
values, a larger contribution comes from high redshift, where 
the space density of clusters is higher in the ACDM model 
(beyond z ~ 0.4 - see Fig. [l]), causing the ACDM curve to 
rise more steeply than the rCDM curve. 

As noted in Bartelmann (2000), the number of sources 
above a given Y is sensitive to the adopted value of the 
baryon fraction. To illustrate this, we have also plotted the 
rCDM counts with a lower baryon fraction, which is in con- 
flict with cluster measurements but consistent with nucle- 
osynthesis determinations. As stated before, using fib = 
0.019/1" 2 (Buries & Tytler 1998) and taking h = 0.5 im- 
plies /b = 0.076. Since Y scales linearly with fb, the factor 
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Figure 3. The number of haloes above a given Y200 (defined 
in the text) for the two surveys. As well as the values obtained 
assuming the cluster baryon fraction, we show the result for 
the tCDM survey with a baryon fraction around a factor of 
two smaller, consistent with standard nucleosynthesis. The ver- 
tical dashed line corresponds to the sensitivity limit of Planck 
(Yu m ~3x 10" 4 arcmin 2 ). 

of ~ 2 difference between the baryon fractions scales the 
counts horizontally by the same factor. However, from Fig. ^ 
the cumulative number of sources scales approximately pro- 
portional to i^oo'' 3 f° r the tCDM case. The change in the 
predicted counts when approximately doubling /b is an in- 
crease by a factor of ~ 2.7. Hence, the degeneracy between 
/b and Q,q implies that the projected source counts is not a 
good constraint on one without accurate prior knowledge of 
the other. This is especially true because the counts are also 
expected to be quite sensitive to the value of as. 



tCDM ACDM 




100 1000 100 1000 
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Figure 4. Source counts above a given absolute flux, in mJy, 
for the rCDM cosmology (left column) and ACDM cosmology 
(right column). The rows from top to bottom specify results for 
central frequencies of 143, 217 and 353 GHz respectively. Within 
each panel, results are plotted for the NFW model (solid lines 
correspond to total absolute fluxes and dotted lines correspond 
to the kinetic fluxes) and /3 model (long-dashed and short-dashed 
lines correspond to total and kinetic fluxes). 



3.2 Source counts 

The number of sources per square degree above a given ab- 
solute flux, in mJy, is plotted in Fig. H for two models: one 
assuming an NFW density profile and one assuming the /3 
model density profile (we will subsequently refer to these 
models as "NFW" and "/3"). The left column illustrates re- 
sults for the rCDM models and the right column for ACDM 
models. From top to bottom, results are for central fre- 
quencies of 143, 217 and 353 GHz respectively, assuming 
a beam FWHM of 8 arcminutes for the 143 GHz results and 
5 arcminutes for the 217 and 353 GHz results. Within each 
panel, the solid and dotted curves are for the NFW model, il- 
lustrating counts for the total (thermal plus kinetic) flux and 
the kinetic component respectively. Long and short dashed 
curves are results for total and kinetic fluxes for the /3 model. 

Comparing the total flux counts between the rCDM 
and ACDM results for the same density profile at 143 GHz, 
the same trend is apparent as in Fig. K at low fluxes both 
models predict similar source counts (e.g. the j3 model pre- 
dicts 0.3 — 0.4 sources per square degree above 100 mJy), 



but at higher fluxes the rCDM model predicts more sources 
than the ACDM model. Using the nucleosynthesis determi- 
nation of the baryon fraction for both models would lower 
the rCDM counts relative to the ACDM counts. 

At 217 GHz, the signal is completely dominated by the 
kinetic effect, since at this frequency the thermal effect gives 
a negligible contribution even when accounting for the fi- 
nite width of the detector's response function. The number 
of sources above a given flux in this channel is at least an 
order of magnitude lower than in the other two channels 
(which can also be seen for a given frequency, when com- 
paring total to kinetic flux), because the distortion caused 
by the kinetic effect is weaker than the thermal effect. As a 
result, the brightest sources observed at 143 and 353 GHz 
have their fluxes substantially reduced at 217 GHz. Further- 
more, since the difference between the rCDM and ACDM 
models is greatest for the largest Y values, the difference in 
source counts between the two cosmologies at 100 mJy is 
greater than at the other frequencies. 

For the 353 GHz channel, the source counts are higher 
at constant total flux than the other two channels. Compar- 
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ing this to the 143 GHz results (i.e. where the signal is also 
dominated by the thermal effect), two differences cause the 
rise. First, the beam FWHM has decreased from 8 arcmin- 
utes to 5 arcminutes which allows more sources to be re- 
solved. Second, for the constant bandwidth (/S.v/v — 0.37), 
the dispersion of the beam a is proportional to v and so the 
SZ flux is higher at 353 GHz than at 143 GHz. 

At all frequencies and in both cosmologies, the j3 model 
predicts more sources above a given flux than the NFW 
model. Since the (3 model density profile is shallower at large 
radii than the NFW profile, corresponding fluxes in the for- 
mer model are larger for individual objects which conse- 
quently increases the total number of sources above a fixed 
flux. 

3.3 Redshift distributions 

Redshift distributions of SZ sources are given in Fig. |J 
which is split into 4 panels. The panels on the left of the 
figure illustrate results for the rCDM survey while the pan- 
els on the right illustrate results for the ACDM survey. 
Each sample only contains haloes with Y > Yu m . The solid 
lines represent the redshift distributions of haloes assum- 
ing Y — I200, with no prior assumption of how the gas is 
distributed within R200 (equation h4j). Dotted and dashed 
lines correspond to distributions for the NFW and /3 models 
respectively. The top panels assume a beam FWHM of 5 ar- 
cminutes and the bottom panels assume a beam FWHM of 8 
arcminutes. (In practice the limit for detection of SZ clusters 
comes from combining different channels which have differ- 
ent FWHM, and we expect the effective beam profile for SZ 
detection by Planck to lie between these values.) 

We first consider the Y200 redshift distributions, which 
assume no beam convolution and in effect assume all clus- 
ters are unresolved by Planck and contribute signal only 
from their inner regions. The amplitude of the peak in the 
redshift distributions for both cosmological models is over 
an order of magnitude lower than in Fig. [l], where haloes of 
all masses are considered. As is evident from Fig. [], there 
is a significant number of sources with Y < Y\i m , and so 
the deficit in clusters in Fig. |H] is due to the removal of these 
sources, and this confirms the near completeness of the sam- 
ple. The peak redshift values for the Y200 distributions are 
z ~ 0.1 for rCDM and z ~ 0.4 for ACDM respectively, 
which are also considerably lower than when the full range 
of haloes are considered. Since Y200 oc M|qq (1 + z)/d 2 A {z), 
the imposed lower limit in Y200 implies a lower mass limit 
which increases with z, due to the effect of d,A- Hence, this 
reduces the number of clusters by a larger factor at high z 
than at low z, causing the peak of the distribution to move 
to a lower redshift. Note in particular that the height of the 
ACDM redshift distribution is now lower than in the rCDM 
redshift distribution, by around a factor of 2.5, since the 
peak in the original ACDM sample is at a higher redshift 
than in the rCDM sample. However, at higher redshift (be- 
yond z ~ 0.4) the ACDM model still predicts an excess of 
sources compared to the rCDM model. 

Comparing the Y200 redshift distributions to those pro- 
duced by the NFW model (in the top panels of Fig. |{| as- 
suming a beam FWHM of 5 arcminutes), the NFW model 
predicts fewer clusters per square degree than the Y200 model 
for redshifts below ~ 0.5 for rCDM (~ 0.7 for ACDM) but 



rCDM ACDM 




0.5 1 1.5 0.5 1 1.5 



z z 

Figure 5. Redshift distributions of haloes in the rCDM sur- 
vey (left column) and ACDM survey (right column), with Y > 
3 X 10 — 4 arcmin 2 . Each panel illustrates results for 3 different es- 
timates of Y, labelled as l^ODi NFW & (see text for details). 
Results for the latter two cases were generated assuming a Gaus- 
sian beam function with a FWHM of 5 arcminutes (top panels) 
and 8 arcminutes (bottom panels). 



more clusters beyond these redshifts. There are several rea- 
sons for the difference. 

Firstly, the size of each cluster for the NFW model 
(i.e. the upper limit of the Y integral) is no longer fixed at 
R200 but is determined by angular radius at which y c — j/b g 
(i.e. #b g ). This allows clusters with sizes 9h s that are less or 
greater than #200 to be included in the sample, providing 
they have Kjfw > Viim- Note also that since y is an integral 
along the line of sight, the NFW model does not calculate y c 
assuming that each halo has a diameter 2.R200 in this direc- 
tion; rather, the limits are selected such that the integrated 
y has converged. 

Secondly, for a halo with constant M200 the value of 
y increases at fixed radius with redshift, because the den- 
sity profile has a higher normalization at higher redshift. 
This is due to the definition of each halo being a region of 
large density contrast (i.e. on average 200 relative to criti- 
cal). As redshift increases, the critical density increases as 
fi _1 (;z) (1 + z) 3 , causing an increase in the average physi- 
cal density of the halo. This causes the radius at which the 
integrated y reaches the background to grow with redshift, 
increasing the value of Y. 

Thirdly, when the y profile is convolved with the beam 
function, it is flattened out to a radius approximately equiva- 
lent to the angular size of the beam. As redshift is increased, 
the apparent angular size of the cluster decreases and so the 
convolution has a significant effect on the profile out to large 
physical radii. 

Hence, the deficit of clusters in the NFW model at low 
z, when compared to the Y200 model, is due to low-mass ob- 
jects which are not resolvable against the background when 
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convolved with the beam (Ynfw < Yu m ). However, at higher 
redshifts the net effect is that the profiles become higher rel- 
ative to the background (even though the convolution has 
a greater effect), leading to a larger angular size and conse- 
quently, for a halo of constant M200, increasing Ynfw rela- 
tive to Y200- 

The same argument can be applied to the f3 model, 
albeit predicting more clusters than the NFW model for ap- 
proximately z > 0.2. The (3 model density profile decreases 
more weakly at large radii than the NFW profile: pp ~ r -2,5 



for our choice of /3 = 5/6, whereas pnfw 



Hence, 



Yp > Ynfw f° r an but the very lowest redshifts where con- 
volution effects dominate. 

In the bottom panels of Fig. |H| results are displayed for 
a beam FWHM of 8 arcminutes; remember that the effective 
resolution of Planck for SZ will lie between 5 and 8 arcmin- 
utes. The number of clusters in both the NFW and f3 mod- 
els decreases when compared to the previous results (by less 
than a factor of 2), although the value of the peak redshift 
is not significantly affected. The larger beam FWHM causes 
the convolution to have an effect out to larger radii, which 
reduces the number of objects able to be resolved against 
the background. 

3.4 Cluster sizes 

While most clusters are effectively point sources at the 
Planck resolution, it is expected to resolve the nearer ones. 
In Fig. we plot the number of clusters, as a function of 
their angular size 8 3 . 

The Y200 model shows the angular sizes corresponding 
to the clusters' 7?200 values, #200 = R20o/<1a{z). The num- 
ber of sources peaks at around 5-6 arcminutes for the rCDM 
case and around 3-4 arcminutes for the ACDM case. Both 
distributions are skewed such that there are more objects 
with sizes greater than the peak size than below it. Broadly 
speaking, the objects with low 8„ values are at high red- 
shift (the minimum resolved angular size of a source with 
Y200 > Yim strongly declines with redshift until z ~ 0.5, af- 
ter which it flattens off), and the objects with large 9 S values 
are at low redshift. The rCDM cosmology predicts similar 
numbers of sources to the ACDM cosmology at the peak of 
the distributions but the former predicts more with high 8 S 
values, due to the excess of sources at low redshift relative 
to the ACDM case. 

However the -R200 radii do not give a good estimate of 
the angular size detectable by Planck, for which we need 
to mode l the cluster profiles using 8 S as described in Sec- 
tion 2.6. We consider beam convolutions of 5 and 8 arcmin- 
utes, as shown by the dotted vertical lines; recall that the 
SZ source catalogue will arise from combining the different 
frequency channels and so the effective Planck resolution in 
SZ is expected to lie between these values. For the NFW 
and j3 models the distributions are much narrower, and the 
peak number of sources higher, than for the Y200 model. The 
sharp cut-off at low 8 S is due to the size at which the con- 
volved high-redshift sources no longer contain enough flux 
to be detected above the sensitivity limit. For a 5 arcminute 
beam (top panels) this size is larger for the f3 model than 
for the NFW model, but for an 8 arcminute beam the two 
are comparable. Thus, the smallest resolved size is sensitive 
to the inner shape of the convolved profiles, determined by 




10 
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Figure 6. Angular size distributions of clusters in the rCDM 
survey (left column) and ACDM survey (right column), with Y > 
3 X 10 _4 arcmin 2 . Solid curves depict distributions assuming Y200 
values, dotted curves are for the NFW model and dashed values 
are for the (3 model. The dotted vertical lines illustrate the value 
of the FWHM, which is 5 arcminutes in the top panels and 8 
arcminutes in the bottom panels. 



both the intrinsic y profile (which is different for the NFW 
and j3 cases) and the beam size. 

The position and height of the peaks in the NFW and (3 
model size distributions indicate that the majority of objects 
will be unresolved by Planck. The largest number of resolved 
sources is found by assuming the rCDM cosmology and a 5 
arcminute beam; even then, the NFW model predicts that 
less than 1 per cent of the sources will be resolved, although 
this rises by almost an order of magnitude for the (3 model. 
This difference is caused by the more extended distribution 
of gas in the /3 model, increasing the number of detected 
sources above the background. We conclude that without 
more detailed information on the cluster profiles it is difficult 
to accurately predict the fraction of sources which will be 
resolvable. 



4 CONCLUSIONS 

In this paper, we have used the 'lightcones' generated from 
the largest dark matter simulations to date, the Hubble 
Volume simulations, in order to predict the distribution of 
sources that will be observable, via the Sunyaev-Zel'dovich 
(SZ) effect, by the forthcoming Planck Surveyor mission. We 
have used the simulations of a critical-density Universe with 
a modified initial CDM power spectrum (rCDM), and a flat 
low-density Universe (ACDM) to obtain masses, redshifts 
and peculiar velocities of dark matter haloes from z — to 
z ~ 4. By making simple assumptions on the distribution of 
ionized gas within these haloes, we have been able to pre- 
dict the expected distribution of SZ sources as a function of 
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redshift, angular size and flux in 3 of the frequency chan- 
nels relevant to the Planck HFI detector. Our calculations 
involve simple modelling of the detector response function, 
beam convolution and the presence of an unresolved back- 
ground. Our main conclusions are as follows. 

The number of SZ sources above a given Y rises at a 
higher rate (for decreasing flux limit) in the ACDM cosmol- 
ogy than for the rCDM cosmology, due to differing redshift 
distributions and the angular diameter distance-redshift re- 
lation. The ACDM model produces more objects at higher 
redshift than does the rCDM model, as well as the angular 
diameter distance being larger for a given redshift. How- 
ever, the absolute number of sources above any given Y is 
not itse lf a good discriminator between differing cosmolog — 



soil 

ick 



ical moHels The normalization is sensitive to the choice of 
baryon fraction (/b), since Y oc /b and N(> Y) falls off 
rapidly with Y, and a significant dependence on ag is also 
expected. 

As expected, the source counts at 217 GHz are domi- 
nated by the kinetic effect, whereas at 143 and 353 GHz the 
counts are dominated by the thermal effect. At 217 GHz, 
the predicted number of sources expected to be detected by 
Planck above 100 mjy is around 0.04 per square degree for 
the rCDM model (with cluster baryon fraction) , lowering by 
about a factor of 5 for the ACDM model. At 143 GHz (as- 
suming a beam F WHM of 8 arcminutes) , this rises to around 
0.3 per square degree for both models, and at 353 GHz the 
narrower beam size implies a value of around 1 source per 
square degree with Y > 3 x 10~ 4 arcmin 2 (appropriate for 
Planck, Bartelmann 2000). Differences in the source counts 
b etween NFW and 8 models are within a factor of 2. 



Ths redshift distribution of SZ sources peaks at z ~ 0.2 



for the rCDM model and z ~ 0.3-0.4 for the ACDM model. 
At redshifts below the peak values, the expected number of 
clusters is insensitive to the adopted gas distribution. How- 
ever, for redshifts above the peak, a more extended gas (pres- 



sure) distribution gives a larger iiumbui uf Semites. Fui uui 
adopte d NFW and (3 donaity profiled, the latter 



v- 



distribution that falls off more slowly with radius, since the 
asymptotic slope of our adopted 3 density profile is r -2 ' 5 , 
whereas the outer slope for the NFW profile is r~ 3 . Increas- 
ing the FWHM of the beam from 5 to 8 arcminutes decreases 
the number of objects. 

Angular sizes of objects also depend on the gas distri- 
butions. For small angular sizes, both 8 and NFW models 
predict a sharp cutoff which is due to the intrinsic shape of 
the convolved y profiles. The vast majority of objects have 
sizes below the beam FWHM and so are not resolved by 
Planck. For the larger angular sizes, the 8 model predicts 
more sources, again due to the more extended distribution 
of gas in the outer regions of the clusters. Both rCDM and 
ACDM models predict similar peak numbers of clusters but 
the rCDM model predicts more resolved clusters due to the 
lower angular diameter distance. 
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